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Abstract. We apply statistical physics to study the task of resource allocation in 
random sparse networks with limited bandwidths for the transportation of resources 
along the links. Recursive relations from the Bethe approximation are converted into 
useful algorithms. Bottlenecks emerge when the bandwidths are small, causing an 
increase in the fraction of idle links. For a given total bandwidth per node, the efficiency 
of allocation increases with the network connectivity. In the high connectivity limit, 
we find a phase transition at a critical bandwidth, above which clusters of balanced 
nodes appear, characterised by a profile of homogenized resource allocation similar to 
the Maxwell's construction. 
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1. Introduction 

Analytical techniques developed in the statistical physics of spin models have been 
widely employed in the analysis of systems in a wide variety of fields, such as neural 
networks [HE], econophysical models [3], and error-correcting codes [21 Hj- Analogy is 
drawn between the interactions and dynamics commonly present in the spin models and 
other systems. Consequently, notions such as cost functions and noises can find their 
corresponding physical quantities in spin models. This parallelism contributes to the 
success of statistical mechanics in various fields. 

Recently, a statistical physics perspective was successfully applied to the problem of 
resource allocation on sparse random networks [SI IB]- Resource allocation is a well known 
network problem in the areas of computer science and operations management [7J [8]. 
It is relevant to applications such as load balancing in computer networks, reducing 
Internet traffic congestion, and streamlining network flow of commodities [91 [10]. In a 
typical setup, each node of the network has its own demand or supply of resources, and 
the task is to transport the resources through the links to satisfy the demands while a 
global transportation cost function is optimized. 

The work in [SI E] can be extended to consider the effects of bandwidths of the 
transportation links. In communication networks, connections usually have assigned 
bandwidths. Bandwidths limit the currents flowing in the links or, in equivalent 
models, the interaction strengths defined on the links connecting the site variables. The 
significance of finite bandwidths was recently recognized in several similar problems. For 
example, transportation problems with global constraints on interaction strengths were 
studied, through a model of transportation network with limited total conductance |HJ, 
and another model of transportation network with limited total pipe volume or surface 
area which increase the resistance of pipes [12]. Among these studies, constraints on 
interaction strengths were implemented as limits on conductance. 

The purpose of this paper is twofold. First, we demonstrate the close relation 
between statistical mechanics and distributed algorithms. Such close relations have 
facilitated statistical mechanics to provide insights to a number of principled algorithms 
recently, as illustrated by the relation between Bethe approximation and error-correcting 
codes [21 H31 [H], and probabilistic inference [IS], as well as the replica symmetry- 
breaking ansatz and the survey propagation algorithm in i^-satisfiability problems [16J. 
These contributions were made in problems with discrete variables. In this paper, we 
extend the Bethe approximation to the resource allocation involving continuous variable 
with finite bandwidths. Since bandwidths limit the currents along the links, both the 
message-passing and price iteration algorithms proposed in [SI E] have to be modified, 
which enable us to find the optimal solutions without the need of a global optimizer. 

The second purpose is to study the behavior of the optimized networks when 
the connectivity and bandwidth changes, using insights generated from the recursion 
relations of the free energy. We observe a number of interesting physical phenomena 
in networks with finite bandwidths. For example, there is the emergence of bottlenecks 
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when the bandwidth decreases, resulting in the shrinking of the fraction of unsaturated 
links. We also identify the correlations between the capacities of the nodes and their 
roles as sources, sinks and relays. We find that resources are more efficiently distributed 
with increasing connectivity. Scaling relations are found in the distributions of currents 
and chemical potentials. In the high connectivity limit, we find a phase transition at a 
critical bandwidth, above which clusters of balanced nodes appear, characterised by a 
profile of homogenized resource distribution reminiscent of the Maxwell's construction 
in thermodynamics. When adapted to scale-free networks, changes in this profile enable 
us to identify the enhanced homogeneity of resources brought by the presence of hubs 
to nodes with low connectivity. 

The paper is organized as follows. In Section we introduce the general model, 
and the analysis is presented in Section [3j We demonstrate the conversion of the derived 
recursive equations into algorithms in Section HI In Section El we compare the analytic 
solutions with numerical simulations in networks with low connectivity, and report on 
the bottleneck effect. We will examine the limit of high connectivity in Section El We 
conclude the paper and point out the potential applications of the model in Section 

2. The model 

We address the problem of resource allocation on random sparse networks of nodes and 
links. Each node of the network has its own supply or demand of resources, and the task 
is to transport the resources through the links to satisfy the demands while optimizing 
a cost function of transportation costs. The amount of resources transported through 
a link is limited by the bandwidth, which leads to extra constraints in the resource 
allocation problem. 

Bandwidth limitations require us to generalize the problem in the following way. 
In the case of infinite bandwidths, currents along the links can be set arbitrarily large 
to satisfy the demands of all nodes. Hence hard constraints on the satisfiability of the 
node demands can be realized, provided that the networkwide supply of resource is 
greater than their resource demand. On the other hand, in the present model of finite 
bandwidth, a node with resource demand can still experience shortage even though its 
neighbors have adequate supply of resources, since the provision of resources can be 
limited by the bandwidths of the links. Hence, the hard constraints on node satisfaction 
of the infinite-bandwidth model is replaced by soft constraints, and the energy function 
is generalized to include the cost of unsatisfaction, or shortage of resources. 

Specifically, we consider a network with N nodes, labelled i — 1, . . . , N. Each node i 
is randomly connected to c other nodes. The connectivity matrix is given by A%j = 1,0 
for connected and unconnected node pairs respectively. Except for the discussion on 
scale-free networks in Section [6] B4, we focus on sparse networks, namely, those of 
intensive connectivity c~ O(l) < iV. 

Each node % has a capacity A$ randomly drawn from a distribution p(Aj). Positive 
and negative values of Aj correspond to supply and demand of resources respectively. 
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The task of resource allocation involves transporting resources between nodes such that 
the demands of the nodes can be satisfied to the largest extent. Hence we assign 
Uij = —i/ji to be the current drawn from node j to i, aiming at reducing the shortage £j 
of node i defined by 



The magnitudes of the currents are bounded by the bandwidth W, i.e., \yij\ < W . For 
simplicity and clarity, we deal with the case of homogeneous connectivity and bandwidth 
in the paper, but the analyses and algorithms can be trivially generalized to handle real 
networks having heterogeneous connectivity and bandwidth. 

To minimize the shortage of resources after their allocation, we include in the total 
cost both the shortage cost and the transportation cost. Hence, the general cost function 
of the system can be written as 



The summation (ij) corresponds to summation over all node pairs, and Aj is a quenched 
variable defined on node i. 

In the present model of resource allocation, the first and second terms correspond 
to the transportation and shortage costs respectively. The parameter R corresponds to 
the resistance on the currents, and A, is the capacity of node i. The transportation 
cost <j>(yij) can be a general even function of y^-. (Our model can be generalized 
to consider transportation costs 4> that are odd functions of yij. In such cases, one 
simply has to introduce extra quenched variables of the links representing the favored 
directions of the currents.) In this paper, we consider and if) to be concave functions 
of their arguments, that is, <f/(y) and ip'(£) are non- decreasing functions. Specifically, 
we have the quadratic transportation cost <f>(y) = y 2 /2, and the quadratic shortage cost 
if>(Ai, {yij\Aij = 1}) = £j 2 /2. As considered previously, other concave and continuous 
functions of 0, such as the anharmonic function in [Bj, result in similar network behavior. 

The model can also be extended to probabilistic inference on graphical models 
[17] . In this context yij = yji may represent the coupling between observables in nodes 
j and i, 4>(yij) may correspond to the logarithm of the prior distribution of y^, and 
if)(Ai, {yij\Aij = 1}) the logarithm of the likelihood of the observables Aj. 

3. Analysis 

The analysis of the model is made convenient by the introduction of the variables It 
can be written as the minimization of Eq. (j2J) in the space of y^ and subject to the 
constraints 




(1) 




(2) 



Ai + ^2 AijVij + ti > 



(3) 
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and the constraints on the bandwidths of the links 
W- Vij >0, 

W + Vij >0. (4) 

We consider the dual of the original optimization problem. Introducing Lagrange 
multipliers to the above inequality constraints, the function to be minimized becomes 

L = [^(&) + ^ ( Ai + Ai i yi i + &) + 

i (ij) 

+ £ A, [mva) + rij(W - Vlj ) + % (W + Vij )] , (5) 

(ij) 

with the Kuhn- Tucker condition 

Hi (Ai + Y AijVij + &) = °» 
(ij) 

a& = 0, 

7S(W-y y )=0, 

l- J (W + y tJ ) = 0, (6) 

and the constraints /ij < 0, oti < 0, 7^ < and 7^ < 0. Optimizing L with respect to 
yij, one obtains 

Vij = y (n - AO ( 7 ) 



y (x) = max< — W 7 , min 



x 



(8) 



The Lagrange multiplier /^ is referred to as the chemical potential of node i, and <fi' is 
the derivative of 4> with respect to its argument. The function Y(/u,j — /ij) relates the 
potential difference between nodes % and j to the current driven from node j to i. For 
the quadratic cost, it consists of a linear segment between fij — Hi = ±WR reminiscent 
of Ohm's law in electric circuits. Beyond this range, y is bounded above and below by 
±W respectively. Thus, obtaining the optimized configuration of currents among 
the nodes is equivalent to finding the corresponding set of chemical potentials /ij, from 
which the optimized y^'s are then derived from Y(fij — /ij). This implies that we can 
consider the original optimization problem in the space of chemical potentials. 

The network behavior depends on whether ip'(0) is zero or not; ip'(0) can be 
considered as the friction necessary for the chemical potentials to overcome when the 
shortage changes from zero to nonzero. For the frictionless case of ip'(0) = 0, to which 
the quadratic ip considered in this paper belongs, two types of nodes can be identified. 
Those nodes % with excess resources (Aj + £V . AijVij > 0) are characterized by /ij = 
and & = 0. Those nodes % with resource shortage are characterized by /ij < and & > 0. 

For the friction case of tp'(0) > 0, there is a third type of nodes % with fully utilized 
resources (Aj + X^^i^i = 0) an d characterized by /ij < and £j = 0. 

Three types of links can be identified. Those links (ij) with jy^l = W are referred 
to as the saturated links. Those with < jy^l < W and \yij\ = are referred to as 
unsaturated and idle links respectively. 
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We introduce the free energy at a temperature T = (3 l , 

F = -T\nZ, (9) 
where Z is the partition function 

f w 

Z = IJ / d y*3 exp ~^ Aii^yii) - I 3 ^( Ai ' Ofal-^i = !}) • 
to) 7_w/ («) 

(10) 

The statistical mechanical analysis of the free energy can be carried out using the replica 
method or the Bethe approximation, both yielding the same results [6]. The replica 
approach for networks with finite bandwidths is a direct generalization of the case 
with infinite bandwidth in [6]. Here, we describe the Bethe approach, whose physical 
interpretation is more transparent. 

In large sparse networks, the probability of finding loops of finite lengths on the 
network is low, and the local environment of a node resembles a tree. In the Bethe 
approximation, a node is connected to c branches of the tree, and the correlations 
among the branches are neglected. In each branch, nodes are arranged in generations, 
A node is connected to an ancestor node of the previous generation, and another c — 1 
descendent nodes of the next generation. 

We consider the vertex V(T) of a tree T. We let F(y\T) be the free energy of the 
tree when a current y is drawn from the vertex by its ancestor node. One can express 
F(y\T) in terms of the free energies Fly^T^) of its descendents k — 1, . . . , c — 1, 

F(y\T) = -TlJ JJ / dy k exp -/3^F(y fc |T fc ) 
L fc=1 \J-w J jfe — i 

c— 1 c— 1 ^ 

k=l k=l ' 

where represents the tree terminated at the k th descendent of the vertex, and Av(t) 
is the capacity of V(T). We then consider the free energy as the sum of two parts, 

F(y\T)=N T F ;iV + F v (y\T), (12) 

where Nt is the number of nodes in the tree T, and F av is the vertex free energy per 
node. Fv(y\T) is referred to as the vertex free energy. Note that when a vertex is added 
to a tree, there is a change in the free energy due to the added vertex. In the language of 
the cavity method [18] , the vertex free energies are equivalent to the cavity fields, since 
they describe the state of the system when the ancestor node is absent. From Eq. ( fl2l) . 
we obtain a recursion relation for the vertex free energy. In the zero temperature limit, 
this relation reduces to 

rC-l 



F v (y\T) = min 

{y k \\yk\<w} 



J2[Myk\Tk) + R(j>(y k ) 
k=i ^ 

+ ^max(-A y (T) -J^yk + y, 0)\ 



(13) 
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Note that F av has to be subtracted from the above relation, since the number of nodes 
increases by 1 when a new vertex is added to the trees. 

The average free energy is obtained by considering the average change in the free 
energy when a vertex is connected to its c neighbors, resulting in 



where Ay is the capacity of the vertex V fed by c trees Ti,...,T c , and (. . .)a 
represents the average over the distribution p(A). The shortage distribution and current 
distribution can be derived accordingly [6]. 

4. The Message-Passing Algorithm 

Due to the local nature of the recursion relation (ITBl . the optimization problem can be 
solved by message-passing approaches, which have been successful in problems such as 
error-correcting cods [13] and probabilistic inference [15]. However, in contrast to other 
message-passing algorithms which pass conditional probability estimates of discrete 
variables to neighboring nodes, the messages in the present context are more complex, 
since they are free energy functions Fy(y\T) of the continuous variable y. Inspired 
by the success of replacing the function messages by their first and second derivatives 
in [6], we follow the same route and simplify the messages devised for this problem. 
Note that the validity of this simplification is not obvious in this problem, owing to 
the non-quadratic nature of the energy function as a result of the finite bandwidths. 
These additional constraints increase the complexity of the problem and complicate the 
messages. 

We form two-parameter messages using the first and second derivatives of the vertex 
free energies. Let (Aij,Bij) = (dF v {yij\Tj)/dyij,d 2 F v {yij\Tj)/dyfj) be the messages 
passed from node j to its ancestor node i, based on the messages received from its 
descendents in the tree Tj. The recursion relation of the message can be obtained by 
minimizing the vertex free energy in the space of the current adjustments drawn 
from the descendents. We minimize 





(14) 



F ij = Yl A i k A i k£ i k + o B i k£ % + R^jktjk + ^"jk^jk + (15) 



subject to the constraints 




(16) 



> 0, 

together with the constraints on bandwidths, 
\Vjk + ejfcl < W. 



(17) 
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(p'- k and 4>"- k represent the first and second derivatives of (f>(y) at y = y^ respectively. 
We introduce Lagrange multiplier for constraints ffl6l) . After optimizing the energy 
function of node j, the messages from node j to i are given by 



Aij 



f'iji 





(18) 



for hrMQ) > 0, 



Ajk{R<t>"jk + Bjk) 1 



for -V»'(0) < /^(0) < 0, 
(19) 



x6 



R(j)' jk + A ifc + 



R( t>"k + 



for ^(O) < -iP'(0), 



where 



9ij(x) = [ip' o %](x) + x, 





^(0) 



for ^(O) > 0, 

for -V»'(0) < /^(O) < 0, 

for h£(0) < -^(0), 



(20) 



(21) 



and hij(x) 



is defined by 



Vij 



Aj — max 

k^i 



—W, min 



R<P'+A jk + xi 



(22) 



These variables can be interpreted as the cavity variables of node j in the absence of i, 
when a current yij is drawn from node j. p^j is the cavity chemical potential of node j. 
hij(x) is the cavity shortage of resource at node j when /z^- takes the value x. if)' o hy{x) 
is then the corresponding dissatisfaction cost per unit resource of node j. The condition 
that if)' o hij(x) + x = in Eq. (|2"T|) is thus to require the chemical potential to be set at 
the minus of the dissatisfaction cost per unit resource shortage. 

For the frictionless case of ip'(0) = 0, if)'(h) changes continuously when h varies 
from zero to negative. Hence when y^ changes, if)' o hij changes continuously, resulting 
in a continuous change in the chemical potential //^ as well as the first derivatives of the 
vertex free energy function A^. Hence, for the quadratic load balancing task, defined by 
<f>(y) = y 2 /2, the vertex free energies in the recursion relation Eq. ffl3l) will be piecewise 
quadratic with continuous slopes, with respect to continuous changes of currents y^, 
implying the consistency of the proposed optimization algorithm. This validities the 
replacement of the message functions by the two-parameter messages. 
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Since the messages are simplified to be the first two derivatives of the vertex free 
energies, it is essential for the nodes to determine the working points at which the 
derivatives are taken. That is, each node needs to estimate the current yij drawn by 
its ancestor. This determination of the working point is achieved by passing additional 
information-provision messages among the nodes. Here, we describe the method of 
backward information-provision messages. This is in contrast to conventional message- 
passing algorithms, in which messages are passed in the forward direction only. An 
alternative method of forward information-provision messages can also be formulated 
following [5]. 

After the minimization of the vertex free energy at a node j, forward messages 
(Aij,Bij) are sent forward from node j to ancestor node i. Optimal currents are 
computed and sent backward from node j to the descendent nodes k ^ i. These 
backward messages serve as a key in information provision to descendents, so that the 
derivatives in the subsequent messages are to be taken at the updated working points. 
Minimizing the free energy ( 1T51) with respect to yjk, the backward message is found to 
be 

y jk - maxj-W 7 , mm [W, y jk - ^ + ^ j j (23) 

When implemented on real networks, the nodes are not divided into fixed 
generations, and there are no fixed ancestors or descendents among the neighbors of 
a particular node. Individual nodes are randomly chosen and updated, by randomly 
setting one of its neighbors to be a temporary ancestor. This results in independent 
updates of the currents y^ and y^ in the opposite directions of the same link. To 
achieve more efficient convergence of the algorithm, we can set y^ = —yji during the 
updates. For the cost function used in this work, we find that this procedure speeds up 
the balance of demands from the two nodes connected by the link. 

An alternative distributed algorithm can be obtained by iterating the chemical 
potentials of the nodes. The optimal currents are given by Eq. (J7J) in terms of the 
chemical potentials fa which, from Eqs. (0Q) and ([6]), are related to their neighbors via 

{0 for h;\0) > 0, 

h;\0) for -f(0) < h-\0) < 0, (24) 
gr\0) for h~\0) < -^(0), 
where hi(x) and gi(x) are given by 

hi(x) = -Aj - } j A ij Y(fj lj - x), 
j 

gi{x) = if) 1 o hi(x) + x, (25) 

with function Y again given Eq. (jSJ). This provides a simple local iteration method 
for the optimization problem in which the optimal currents can be evaluated from the 
potential differences of neighboring nodes. Simulation results in Section [5] show that 
the chemical potential iterations have an excellent agreement with the message-passing 
algorithm in Eqs. AH)- (EE). 
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We may interpret this algorithm as a price iteration scheme by noting that the 
Lagrangian in Eq. (jSJ) can be written as 



Therefore the problem can be decomposed into independent optimization problems, each 
for a current on a link, /ij is the storage price per unit resource at node i, and each 
problem involves balancing the transportation cost on the link, and the storage cost at 
node i less that at node j. From Eqs. (I!MI) and (|25|) . we see that when node i is short 
of resources, /ij = —ip' o h{pi). That is, the storage price at a node is equal to minus 
the dissatisfaction energy of that node. This means that the algorithm gives bonuses 
to the link controls instead of charging them, so as to encourage them to station their 
resources at the nodes with resource shortage. The amount of bonus per unit resource 
compensates exactly the dissatisfaction energy per unit resource. This corresponds to a 
pricing scheme for the individual links to optimize, which simultaneously optimizes the 
global performance. The concept is based on the same consideration as those used in 
distributed adaptive routing algorithms in communication networks [19]. 

5. Results and Physical Implications 

We examine the properties of the optimized networks by first solving the theoretical 
recursive equation (TT3]) numerically to obtain various quantities of interest, including the 
energy, the fraction of idle links and saturated links. The results are then compared with 
the simulation results obtained by the message-passing and price iteration algorithms. 
All experiments assume a quadratic shortage cost = £ 2 /2 for shortage £ given in 
Eq. (jHJ), and a quadratic transportation cost <p(y) = y 2 /2, and the capacity distribution 
p(A) a Gaussian with mean (A) and variance 1. 

To solve numerically the recursive equation ( 1T31) . we have discretized the vertex 
free energy functions Fy(y\T) into a vector, whose i th component is the value of the 
function corresponding to the current At each generation of the calculation, a node j 
is randomly chosen and its capacity Aj is drawn from the distribution p{Af). The node 
is randomly connected to c — 1 nodes of the previous generation with vertex free energies 
Fv(yk\Tk)- Then the vertex free energy Fy(y\T) is found by minimizing the right hand 
side of Eq. (fl3j) using standard numerical methods for every discretized components 
of y, subject to the constraints of \yk\ < W. Thus a new generation of vertex energy 
functions is generated and the process continue in a recursive manner. 




(26) 



where 




(27) 



5.0.1. Effects of Bandwidths Figure [T] shows the average energy as a function of the 
bandwidth W, obtained by, respectively, the message-passing and the price iteration 
algorithm simulations, and the recursive equation ( [131) . Simulation results collapse 
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Figure 1. Simulation results of average energy per node as a function of bandwidth 
W for N = 1000, c = 3, R = 0.1 and (A) = 0, -0.1, -0.2 with 200 samples, as 
compared to the numerical solution of Eq. (|13[) . x's correspond to the energy in the 
zero bandwidth limit as obtained from Eq. (|28p. Inset: The fraction of idle, unsaturated 
and saturated links as a function of bandwidth for (A) =0. x corresponds to /idle in 
the zero bandwidth limit as obtained from Eq. ([29]) . 



almost perfectly with the recursive equation. We find that the average energy increases 
with decreasing bandwidths. Physically, decreasing bandwidth corresponds to the 
increasing limitations in the allocation of resources. Thus, it is natural to have higher 
energy. In the limit of zero bandwidth with small finite connectivity, the solution 
is equivalent to the initially unoptimized condition. As the bandwidth vanishes, no 
resources are transported among the network, and the average network energy is given 
by 



E = N 



A 2 

— p(A)dA 



(2f 



In the limit of zero bandwidth, a link can only be idle or saturated, with no 
intermediate unsaturated state. A link is idle only if the two connected nodes are 
satisfied already, with their own initial resources. The probability of finding an idle link 
is thus determined by the initial resources of the connected nodes, yielding the fraction 
of idle links as 



idle 



P (A)dA 



1 — erf 



(A) 
y/2 



(29) 



where the last equality is specific to the Gaussian p(A). The fraction of saturated link is 
indeed given by 1 — /idi e - The corresponding limit is shown in the inset of Fig. (TJ which 
is consistent with the simulation results. 



5.0.2. Bottleneck Effect Another interesting physical consequence of decreasing 
bandwidth is also shown in the inset of Fig. [TJ We note that as the bandwidth decreases, 
the fraction of saturated link increases, which is a direct consequence of the attempt to 
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Bandwidth W = 2 Bandwidth W = 1 

Figure 2. An example of bottleneck effect. 

minimize the shortages of nodes fed by the saturated links. A similar agreement would 
lead us to anticipate a decreasing fraction of idle links as the bandwidth decreases, since 
more links should participate in the task of resource allocation. Surprisingly, we notice 
an increasing fraction of idle links as the bandwidth decreases, in contrast with our 
anticipation. 

This is a consequence of the bottleneck effect, as illustrated in Fig. [2j When the 
bandwidth decreases, resources transferred from the secondary neighbors may become 
redundant since resources from nearest neighbors already saturate the link to the 
dissatisfied node, which can therefore be considered as a bottleneck in transportation. 
This makes the resource distribution among the nodes less uniform, since resources are 
less free to be transported through the links. As a result, more nodes are having either 
excess resources or large shortage, leading to a higher network energy. On the other 
hand, the existence of these bottlenecks confirm the efficiency of our algorithms, since 
redundant flows to distant nodes are eliminated. 

The existence of bottlenecks is common in many real networks. Among the most 
common examples are bottlenecks occurring in traffic congestions. Further studies can 
be carried out to minimize the bottleneck effect, such as using heterogeneous bandwidths 
for different links. The highway is an example of enlarging the bandwidth on a link 
with heavy traffic. Similar considerations have been applied to routing of traffic on 
telecommunication networks [T9] . 

5.0.3. Density of Saturated Links and Its Physical Implications To study the roles 
played by the different nodes in resource allocation, we consider the example of networks 
with c = 3 and very negative (A) , so that all resources in the networks have to be utilized 
for optimization. We classify the nodes into four classes, the nodes of each class having 
^sat = to 3 saturated links connected to them. Figure [3] shows the conditional capacity 
distribution of these four classes. 

For nodes having no saturated links, their participation in the optimization task is 
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Figure 3. The distribution of initial resource A given different number of saturated 
links nsat after optimization for N = 1000, (A) = —5, c = 3, R = 0.1 and W = 1 with 
4000 samples. 

relatively inactive, and their capacity distribution is approximately Gaussian with an 
average (A) . For nodes having one saturated link, two peaks are found around A — W 
and A + W in the conditional capacity distribution, respectively corresponding to nodes 
acquiring resources from or providing resources to their neighbors, in order to lower the 
global shortage cost. 

For nodes having two saturated links, we obtain a trimodal distribution. The right 
and the left peaks at around A + 2W and A — 2W correspond to resource donors and 
receptors respectively, which export or import resources, saturating two of the connected 
links. The central peak at A corresponds to nodes with initial resources close to the 
average value. These nodes have little tendency to donate resources to their neighbors 
or receive from them. Thus they are referred as relays which participate in optimization 
by transferring resources from one neighbor to another. For nodes with three saturated 
links, we observe a four-modal distribution. Similarly, the rightmost and the leftmost 
peaks correspond to sources and sinks of resources respectively. The middle-right and 
the middle-left peaks correspond to source-like relays and sink-like relays respectively, 
namely, they serve partially as relays and partially as sources and sinks. All these 
distributions show the strong dependence of the number of saturated links of a node 
on its own initial resources, which identify an inborn optimization role for each node 
during the resource allocation. 

6. The High Connectivity Limit 

The optimal network behavior can be studied analytically in the limit of high 
connectivity. In this limit, the bandwidth of the links plays an important role in the 
scaling laws of the physical quantities of interest. Two cases will be considered in this 
section. In the first case, the bandwidths of the individual links remains constant when 
the connectivity scales up. This increases the total available bandwidth connecting a 
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node, and hence the freedom in resource allocation. 

In the second case, we scale the bandwidth of the individual links as c~ l when the 
connectivity c scales up. Hence, the total available bandwidth connecting the individual 
nodes is conserved, but the burden of transporting resources is divided into smaller 
currents shared by a larger number of neighbors. 

To further simplify the analysis, we set <f>(y) = y 2 /2, ip(£) = £ 2 /2 and R ~ 0(1) in 
the subsequent derivations. 



6.1. High Connectivity with Non-vanishing Bandwidths 

6.1.1. The Optimal Solution In the high connectivity limit, the magnitudes of the 
currents in individual links are reduced, owing to the division of the transportation 
burden among the links. Since the bandwidths remain constant, the fraction of saturated 
links becomes negligible. Hence, Eq. ffT9l) converges to the result 

B tj = e(-fiij)^. (30) 
At the steady state, all current adjustments €jk vanish, reducing Eq. fl23|) to 

Vjh = ^ {Aj - A jk ) . (31) 
Thus, Eq. (1211) for ^ in the message-passing algorithms can be simplified to 

= maxiYl + -0 [ Vij -^j + ^Y. A jk A jk j , j , (32) 

where we have approximated c — 1 by c. This expression of A^ corresponds to the 
shortage on node j after optimization for the quadratic cost = £. We then utilize 
the fact that the messages from the descendents are independent to each other. This 
allows us to express the collective effects of the descendents on a node in terms of the 
statistical properties of the descendents. By virtue of the law of large numbers, it is 
sufficient to consider the mean tjia and the variance cr\ of the messages A^. In Eq. ( |32|) . 
the term i/ij is negligible in the limit of high connectivity. After applying the law of 
large numbers to the term J2 k -AjkAj k , we can write 

' cm a — RAj 
c + R 

Averaging over Aj, which is drawn from a Gaussian distribution of mean (A) and 
variance 1, we obtain a self-consistent expression for tua, 



A i:j = max ( — - \ ) . (33) 



(A) + m A = —j= exp 



1 (cm a 



(A) 



_ ( -_< A)M - 2 _4 (34) 

where H(x) — |[1 — erf(a;/v / 2)]- In the limit of high connectivity, this further implies 
m A = max(-(A),0). (35) 
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This equation relates the average shortage after optimization to the average initial 
resources (A) of the system, implying a uniform resource redistribution in the large 
connectivity limit. Given a negative (A), the equation implies no nodes are keeping 
excess resources relative to the average (A). All excess initial resources are redistributed 
to achieve a global optimized state. We note that the currents scale as (c + R)^ 1 , 
resulting in a much smaller transportation cost compared with the shortage cost when 
the connectivity increases. This allows the re-allocation of initial resources through 
vanishing currents to minimize the shortage cost. The sum of c currents of 0[(c + R)^ 1 ] 
results in a total resource transport of 0(1) transported to or from a node. It provides an 
extremely large freedom in resource allocation. For positive (A), the average shortage 
tua vanishes in the limit high connectivity limit, corresponding to a zero fraction of 
unsatisfied nodes. Hence, increasing connectivity results in better resource allocation of 
the system. 



6.1.2. Network Properties We obtain the current distribution by combining the 
expression for Ay in Eq. fl33l) with the equation for current in Eq. (l3Tj) . such that 



(cm A 
max - 



c + Rl V R 



A,,0 



( cmA a n 
max ( — A j , 



V R 



(36) 



Hence the current distribution is given by 

P(y) = J dAip(Aa) J dA 2 p(A 2 )5^y - 
1 fcrriA K r,\ fcmA 



c + R 



max 



R 



Ai,0 



max 



R 



Ao,0 



For the Gaussian distribution p(A), we obtain 



P(y) = H\ X )KV) + \l ~(c + R)H(X) exp 



--((c + R)y- x ) 2 



(37) 



(38) 



+ c_+Jl H ((c + R)y-2x 



exp 



-~((c + R)V) 



'71 v V2 

where x = cuia/R — (A). This shows that the rescaled distribution P[(c + R)y]/(c + R) 
is independent of c and solely dependent on (A). This confirms that the currents scale 
as (c + -R) -1 as the connectivity increases, as described before. 

We obtain the chemical potential distribution by considering the chemical potential 
of a central node fed by c descendents. Introducing Lagrange multipliers, the 
Lagrangian becomes 

2 



L = ^ + p:[A Q + J^Voj + £o) +^ + J2 [-Y 1 + AojVoj 

3=1 3=1 



+ 0(-^.)^ + XqAW - y 0j ) + XqAW + y 0j 



2c 



(39) 



The currents are given by yoj = (— A j — p)/R, and the chemical potential is given by 



v- = Y^n min ( A o - ^ A v> 1 

R j=i 



(40) 
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Using the statistical properties of A j in the high connectivity limit, the approximate 
expression for fi becomes 



1 • (k aUA nl 
mm Aq U ) . 



R 



(41) 



For Gaussian distribution of the initial resource, we derive the chemical potential 
distribution given by 



1 + - 



cxp , 

2tt 12 



1 



1 + R' fI + X 



Finally, the average energy per node is given by 
(E) = \{^) + f <2/ 2 >, 



e(-//)+if( X )5(//). (42) 



(43) 



where the first term corresponds to the shortage energy per node and the the second 
term corresponds to the transportation energy on the links connected to a node. Using 
Eqs. ( 136]) and (f4T|) . the average energy per node becomes 



(E) 



R 



2(c + R) 



f cmA a n 
max A, 

V R 



Rc 



2(c + R) 2 \ 



( max 



(cm a 
~R~ 



A,0 



For the Gaussian initial resource distribution, this reduces to 



where 



(E) = 

h(x) 
h{x) 



R 



2(c + R) 



h{x) 



Rc 



2(c + R) 2 



hixf 



Dz{\ - z) 



X 

e 2 



^7T 



+ XH{-X); 



X 

e 2 



Dz( X - z) 1 = X-j= + (X 2 + I)H(-X)- 
v 2tt 



(44) 
(45) 

(46) 
(47) 



6.1.3. Efficient Resource Allocation in The Shortage Reigrae In the regime of negative 
(A), apart from the uniform resource redistribution deduced from Eq. fl35l) . the current 
distribution, the chemical potential distribution and the energy also show universal 
characteristics of efficient resource allocation. When (A) is negative, x ^ 1 m the high 
connectivity limit, the current distribution reduces to 



P{y) 



c + R 



7T 



■ exp 



(c + R)Y 



(48) 



Physically, the disappearance of the term 5(y) in the expression of the current 
distribution for negative (A) implies that none of the nodes are holding excess resources 
relative to the average (A), reproducing our conclusion from Eq. fl35l) . 
The chemical potential distribution becomes 



exp 



n 2 



l + ^)(^-(A» 



(49) 
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This shows that the distribution P[(l + c/R){ji — (A))]/(l + c/R) is independent of c 
and depends solely on (A). Thus, the width of the distribution scales as (1 + c/R)' 1 . 
The variable —\i = £q, corresponding to the shortage of a node. This implies that 
the variance of the shortage distribution also scales as (1 + c/R)' 1 , indicating a more 
efficient allocation of resources. The absence of the delta function component at \i — 
implies that none of the nodes are holding excess resources in the negative (A) regime. 

As both the current and chemical potential distributions are Gaussian, the average 
energy per node is simplified to 



Note that in the high connectivity limit, the average energy per node approaches (A) 2 /2. 
This corresponds to the theoretical limit of efficient and uniform resource allocation, in 
which the total shortage is evenly shared among all nodes. 

6.1.4- Comparison with Simulations Simulation results are compared with the 
analytical prediction in the high connectivity limit. Figure H] shows the rescaled current 
distribution P[(c + R)y]/(c + R). For positive (A), the distributions for different 
connectivities collapse almost perfectly after the currents are rescaled by (c + R)' 1 , 
with a very weak dependence on c, approaching the high connectivity limit progressively. 
For negative (A), the distributions also approach the high connectivity limit when the 
connectivity increases, although the dependence on c is larger. 

The distribution has a delta function component at y — and a continuous 
component when (A) is positive. However, when (A) is negative, the delta function 
component disappears and the distribution approaches a Gaussian. As shown in the 
inset of Fig. HI the magnitude of the delta function component, that is, the fraction of 
idle links, vanishes as (A) becomes more negative. Another delta function component is 
observed at \y\ = W when c is small, and dissappears as c increases when the distribution 
approaches the high connectivity limit. 

Figure [5] shows the rescaled chemical potential distributions P[(l + c/i?)(/i — 
(A))]/(l+c/i?) compared with the high connectivity limit. For negative (A) , the rescaled 
distributions are Gaussian-like and approaching the high connectivity limit. Compared 
with the corresponding rescaled current distribution, there is a larger dependence on c. 
The inset of Fig. [5] shows the variances of chemical potentials which approach different 
asymptotic values values as (A) decreases, for different fixed values of connectivity. 
Remarkably, for sufficiently negative (A), the variances become (A) independent. 
This further confirms the convergence of chemical potential distributions to universal 
Gaussian distributions in the negative (A) regime. 

Figure [6] shows the average shortage, as measured by m^, obtained from simulations 
at different connectivities as compared to the theory. In the neighborhood of (A) ps 0, as 
the connectivity increases, tua gradually approaches the high connectivity limit which 
corresponds to most uniform allocation of resources given by Eq. (|35|) . Away from 
(A) « 0, rriA collapses well with the analytical results. 




(50) 
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(c+R)y 



Figure 4. The rescaled current distribution P[(c + R)y]/(c + R) for N = 10000, 
R = 0.1 and W — 1 with 200 samples at different values of (A), as compared to the 
corresponding theoretical large connectivity limits given by Eqs. ([38]) and (I48[) . Inset: 
The fraction of idle links P(y = 0) as a function of (A). 
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Figure 5. The rescaled chemical potential distribution P[(l+c/ R)(fi— (A))]/ (1+c/R) 
for N = 10000, R = 0.1, W = 1 and (A) = -5 with 200 samples, as compared to the 
theoretical large connectivity limit given by Eq. (|49|) . Inset: The variance of chemical 
potential as a function of (A). 



Remarkably, the distributions of currents and chemical potentials, and the average 
shortage m^, all approach their high connectivity limits, for relatively low values 
of c already. We also analyze the finite size effects of the system in approaching 
the asymptotic behavior of the high connectivity limit. The inset of Fig. [6] shows 
the difference between the variance of chemical potentials and its asymptotic value 
(1 + c/R)~ 2 When the connectivity is low, the differences are larger as indicated by the 
power law fit. 

To summarize, our analysis shows that in the high connectivity limit, the optimized 
state of the system is one in which the resources are uniformly allocated, such that 




<A> 



Figure 6. The parameter tua as a function of (A) for N — 10000, R — 0.1 and W = 1 
with 200 samples, as compared to the theoretical large connectivity limit given by 
Eq. (|35p . Inset: The difference <r 2 — (1 + c/R)~ 2 between the variance of the chemical 
potential and its corresponding large connectivity limit, given by Eq. (f4"5|) as a function 
of system size N. The dashed lines shows the power law fits. 



the final resources on every individual node are equal. Simulations with increasing 
connectivities reveal the asymptotic approach to this uniform limit, and the deviation 
from the limit decreases as the connectivity and the system size increase. 



6.2. High Connectivity with Vanishing Bandwidth 

We now consider the case that the bandwidth of individual links scales as W/c when the 
connectivity increases, where W is a constant. Thus the total bandwidth W available 
to an individual node remains a constant. This is applicable to real networks in which 
the allocation of resources is limited by the processing power of individual nodes. 



6.2.1. The Well-defined Chemical Potential Function We start by writing the chemical 
potentials using Eq. (pH|) . 

N 

Hi = min Aj + AijY((j,j - /^), . (51) 
j'=i 

In the high connectivity limit, the interaction of a node with all its connected neighbors 
become self-averaging, making it a function which is singly dependent on its own 
chemical potential, namely, 

N 

AjY{^ - Mi) ~ cMifii). (52) 

i=i 

Physically, the function M(//) corresponds to the average interaction of a node with its 
neighbors when its chemical potential is \l. Thus, we can write Eq. (|5ip as 

fi = min[A + cM(//),0], (53) 
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where /x is now a function of A, and we have 

POO 

M(jn)= / dAp(A)Y(fi(A)- ^ (54) 



where we have written the chemical potential of the neighbors as fi(A), assuming that 
they are well-defined functions of their capacities A. 

To explicitly derive M(/i), we take advantage of the fact that the rescaled 
bandwidth, W/c vanishes in the high connectivity limit, so that the current function 
Y(fij — /ij) is effectively a sign function, which implies that the current on a link is 
always saturated. (This approximation is not fully valid if c is large but finite, and will 
be further refined in subsequent discussions.) Thus, we approximate 

M(m) = — / dAp(A)sgn[//(A) - /4 (55) 

Assuming that //(A) is a monotonic function of A, then we have sgn[/x(A) — fit] = 
sgn(A — Aj), such that 



M{m) = — 

c 



A, 



dAp(A)- / dAp(A) 



(56) 



For Gaussian distribution of capacities, M(//j) is given by 

M W = -ferf(^>). (57) 



Thus, //(A) is explicitly given by 



\i = mm 



(5* 



This equation relates the chemical potential of a node, i.e. the shortage after resource 
allocation, to its initial resource before. Resource allocation through a large number 
of links with a fixed total bandwidth results in a well-defined function relating the two 
quantities. 



6.2.2. Maxwell's Construction Eq. ( 1581) gives a well-defined function /i(A) as long as 
W < a/tt/2. However, a phase transition takes place when W reaches a/tt/2. When 
W > a/tt/2, turning points exists in p(A) as shown in Fig. [7(a), due to the dominant 
effect of the error function component in the equation. p(A) is no longer a monotonic 
function. This creates a thermodynamically unstable scenario, since in the region of 
//(A) with negative slope, nodes with lower capacities have higher chemical potentials 
than their neighbors with higher capacities. This implies that the current flows from 
poorer nodes to richer ones. Mathematically, the non-monotonicity of /x(A) means that 
sgn[/i(A) — fa] and sgn(A — Aj) are no longer necessarily equal, and Eq. fl58l) is no longer 
valid. 

Nevertheless, Eq. (1511) permits another solution of constant p in a range of A. This 
is possible since for any \ii an d Hj m this range, — Hj\ — RW/c implies that the 
link between nodes i and j is unsaturated. If an extensive fraction of the c links of 
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Figure 7. Maxwell's construction on //(A). 



such a node are connected to other nodes in the range, then the freedom of tuning the 
currents in the unsaturated links enables the nodes in the range to have the same level 
of resources after optimization. 

Hence, we propose that the unstable region of //(A) should be replaced by a 
range of constant /i as shown in Fig. 0(b) analogous to Maxwell's construction in 
thermodynamics. Let (A<,yU ) and (A>,/x D ) be the end points of the Maxwell's 
construction as shown in Fig. [7(b). The position of this construction can be determined 
by the conservation of resources. Since the construction is not necessary for positive 
(A), we focus on the case of negative (A). In the high connectivity limit, resources are 
so efficiently allocated that the resources of the rich nodes are maximally allocated to 
the poor nodes as much as the total bandwidth W allows. Let A Q be the capacity of a 
node beyond which nodes can have positive final resources (i.e. zero chemical potential). 
Nodes with A > A Q send out their resources without drawing inward currents from their 
neighbors, and can be regarded as donors. The value of A Q is given by the self-consistent 
equation 

A n = W I dAp(A). (59) 



(60) 



For Gaussian distribution of capacities, 



A, 



W 

Y 



1 + erf 



A, 



-<A) x 

y/2 ) 



By equating the networkwide shortage to minus the initially available resources, we find 



= ^[-A,6(A - A,) - A 6(Aj - A ) 



(61) 



On the right hand side, the first term corresponds to the sharing of shortages among the 

non-donors, while the second term corresponds to the resources donated by the donors 

to reduce the total shortage. In the high connectivity limit, this equation becomes 

"A> / M r A <\ 

dAp(AK- / + / dAp(A)/x(A) 



dAp(A)A-A / dAp(A). 



(62) 



A„ 
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Figure 8. The dependence of the fraction of balanced nodes on the bandwidth W. 



Substituting Eqs. ( 1531) . ( 1561) and ( 1591) . and making use of the symmetry relation 

W / dAp(A) / rfA / p(A , )sgn(A' - A) = (63) 

J oo J oo 

we arrive at the condition 

M> M> 

/i D / dAp(A) = / dAp(A)/i(A), (64) 

</A< </A< 

which implies that the value of p, Q should be chosen such that the areas A and B in 
Fig. [3(b), weighted by the distribution p(A), should be equal - in direct correspondence 
with the Maxwell's construction in thermodynamics. 

For capacity distributions p(A) symmetric with respect to (A), we have 

/x =(A) = i(A < +A > ). (65) 

As a result, the function /^(A) is given by 

(A) for /i< < fi < /i>, 



MA) 



mm 



A - Werf ( A ) , 



(66) 

for otherwise, 



x = (A) + Werff ^ ^ ) . (67) 



where as A< and A> are respectively given by the lesser and greater roots of the equation 

'x - (A)- 

Nodes i with chemical potentials fii = (A) represent clusters of nodes interconnected 
by an extensive fraction of unsaturated links, which provides the freedom to fine tune 
their currents so that the shortages among the nodes are uniform. They will be referred 
to as the balanced nodes. On the other hand, nodes outside the balanced clusters are 
connected by saturated links only. 

The fraction f hal of balanced nodes is given by the equation 

/M = erf(^). (68) 
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Figure 9. (Colour online) The simulation results of fi(A) for TV = 10000, c = 15, 
R = 0.1, (A) = — 1 and W — 3 with 70000 data points, compared with theoretical 
prediction. Inset: The corresponding results for W = 1.2. 

Note that f hal has the same dependence on W for all negative (A). Fig. [8] shows 
that when the total bandwidth W increases, the fraction of balanced nodes increases, 
reflecting the more efficient resource allocation brought by the convenience of increased 
bandwidths. When W becomes very large, a uniform chemical potential of (A) 
networkwide is recovered, converging to the case of non-vanishing bandwidths. 

The chemical potential distribution can be obtained from Eq. fl66l) . Unlike the 
corresponding distribution in the case of non-vanishing bandwidths described in the 
previous subsection, there is no explicit scaling on the connectivity. In addition, it 
is not purely Gaussian. The delta function component 5{n — (A)) corresponds to the 
balanced clusters. The weight of the delta function component is a measure of the 
efficiency of resource allocation. 

We compare the analytical result of /i(A) in Eq. ( 1661) with simulations in Fig.[9j For 
W > data points (A, fi) of individual nodes from network simulations follow the 

analytical result of /x(A), giving an almost perfect overlap of data. The presence of the 
balanced nodes with effectively constant chemical potentials is obvious and essential to 
explain the behavior of the majority of data points from simulations. Outside the region 
of balanced nodes, the data points follow the tails of the function /u(A). For W < 
the analytical /i(A) shows no turning point as shown in the inset of Fig. [9j Despite the 
scattering of data points, they generally follow the trend of the theoretical /i(A). This 
scattering effect may be explained by the use of a finite connectivity in simulations. We 
found that the extent of scattering can be much reduced by increasing the connectivity 
in the simulations. 
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6.2.3. Large but Finite Connectivity Our analysis can be generalized to the case of 
large but finite connectivity, where the approximation in Eq. fl55|) is not fully valid. 
This modifies the chemical potentials of the balanced nodes, for which Eq. ( 1551) has to 
be replaced by 



W 



A< 



M(ji) = — / rfAp(A) - / dAp(A) + / dAp(A) 



oc 



MA) ~ 
R 



(69) 

We introduce an ansatz of a linear relationship between p and A for the balanced nodes, 
namely, 

/j = mA + b. (70) 

After direct substitution of Eq. ( J70|) into M(p) given by Eq. ( 1691) . we get the self- 
consistent equations for m and b, 

R 

m = /, 777T , (71) 



A>-(A) 
V2 



# + c erf ( 



(A). (72) 



Thus, the Maxwell's construction has a non-zero slope when the connectivity is finite. 

We remark that the approximation in Eq. (16"9"|) assumes that the potential differences 
of the balanced nodes lie in the range of 2RW/c, so that their connecting links remain 
unsaturated. Note that the end points of the Maxwell's construction have chemical 
potentials (A) ± RW/c respectively, rendering the approximation in Eq. fl69|) exact at 
one special point, namely, the central point of the Maxwell's construction. Hence, this 
approximation works well in the central region of the Maxwell's construction. However, 
when one approaches the end points of the Maxwell's construction, the balanced nodes 
are also connected to nodes outside the balanced clusters with potential differences less 
than 2RW/c. Hence, we expect to see deviations from the theoretical linear prediction 
near the end points of the Maxwell's construction. 

The chemical potential distribution can be obtained. Compared with corresponding 
distribution in the high connectivity limit, the delta function component of the balanced 
nodes is now smeared into a Gaussian component lying in a strip of width m(A>— A<). It 
implies a lower efficiency in resource allocation due to the finiteness of the connectivity. 

In the simulation data shown in Fig. [TUl the data points of (A, p) from different 
resistance R follow the trend of the corresponding analytical results, both within and 
outside the linear region, with increasing scattering within the linear region as R 
increases. As expected, there are derivations between the analytical and simulational 
results at the two ends of the linear region, with smoothened corners appearing in the 
simulation data, especially in the case of R = 2. We note that when R increases, the 
gradient of the linear region increases, corresponding to a less uniform allocation of 
resources caused by higher transportation costs during optimization. As shown in the 
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Figure 10. (Colour online) Simulation results of (A, y) for N — 2 x 10 5 , W — 3, c = 12 
and (A) = —5 at different values of R, each with 65000 data points, as compared to 
the theoretical predictions. Inset: the corresponding chemical potential distribution 
P(/ii) of the 3 cases. 

inset of Fig. [TDJ the chemical potential distributions follow the trend of the analytical 
results. 

Remarkably, as evident from Eq. (1711) . even with constant available bandwidth W, 
increasing connectivity causes m to decrease, and hence sharpens the chemical potential 
distribution. The narrower distributions correspond to higher efficiency in resource 
allocation. It leads us to realize the potential benefits of increasing connectivity in 
network optimization even for a given constant total bandwidth connecting a node, 
despite a decrease in bandwidth on individual links. 

6.2.4- Resource Allocation in Scale-free Networks Recent studies show that complex 
communication networks have highly hetergeneous structures, and the connectivity 
distribution obeys a power law [20]. These networks are commonly known as scale- 
free networks and are characterized by the presence of hubs. The presence of these 
nodes with very high connectivity can modify the network behavior significantly. 

Figure [TT] shows the simulation results of node i with q = 3 and Cj = 10 in a scale- 
free network, where Cj is the connectivity of node i and is drawn from the distribution 
-P(cj) ~ c~ 3 when the scale-free network is constructed. The data points of (A, /x) follow 
the corresponding analytical results of Eqs. ( 17T1) and ( 1721) . for both sets of nodes with 
Q = 3 and c, = 10. This implies that the previous argument of increasing efficiency 
by increasing connectivity also hold for scale-free networks, as a smaller gradient m is 
found for nodes with higher connectivity. As can be seen from Fig. [TTJ the data points 
are less scattered for nodes with large q, implying a more efficient resource allocation 
for the hubs, in addition to the effect of smaller m. 
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Figure 11. (Colour online) Simulation results of (A,//) for N — 2 x 10 5 , R — 0.1, 
T4^ = 3 and (A) = — 1 as compared with theoretical results, for (a) nodes with Ci = 10 
in scale-free networks with P(c;) ~ c^ 3 , (b) nodes with c$ = 3 in scale-free networks 
and nodes in regular networks with c = 3. Each data set contains 2500 data points. 



More important, nodes with low connectivity benefit from the presence of hubs in 
the networks. To see these benefits, the simulation results of nodes in scale-free networks 
are compared with nodes in regular networks of the same connectivity. As shown in Fig. 
[TlTb). the data points from regular networks are more scattered away from the Maxwell's 
construction, when compared with those from scale-free networks. This shows that the 
presence of hubs increases the efficiency of the entire network, especially for nodes with 
low connectivity. This provides support that scale-free networks are better candidates 
for resource allocation than regular networks. 

7. Conclusion 

We have applied statistical mechanics to a system in which magnitude of the interactions 
(such as currents in resource allocation) between on-site variables (such as the shortages) 
are constrained. This allows us to study an optimization task of resource allocation on a 
sparse network, in which nodes with different capacities are connected by links of finite 
bandwidths. Analogy can be drawn between interaction of spins in magnetic systems 
and that of resources of nodes in networks. The bandwidths serve as constraints on 
the interaction magnitudes and limit the information exchange among the neighbors. 
Despite these additional constraints, we found that the analyses of the unrestricted case 
in [5l E] are applicable after appropriate adaptations, such as allowing for shortages with 
finite penalty and obtaining a non-linear relationship between currents and potential 
differences. 

By adopting suitable cost functions, such as quadratic transportation and shortage 
costs, the model can be applied to the study of realistic networks with constrained 
transport between neighbors. In this paper, we focus on cases in which the shortage 
dominates. By employing the Bethe approximation or equivalently the replica method, 
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recursive relations of the vertex free energies can be derived. 

Analytically, the recursive relations enable us to make theoretical predictions of 
various quantities of interest, including the average energy, the current distribution, and 
the chemical potential distribution. The predictions are confirmed by simulation results. 
In particular, the study reveals interesting effects due to finite bandwidths. When the 
bandwidth decreases, resource allocation is less efficient, and links are more prone to 
saturation. A consequence is the creation of bottlenecks, which refer to the saturation 
of links feeding the poor nodes, rendering the secondary provision of resources from 
their next nearest neighbors redundant. This causes certain links previously assigned 
for secondary transport to become idle, lowering the participation of individual links in 
global optimization and making the resources less uniformly distributed. In the context 
of resource allocation, further studies can be carried out to suppress the bottleneck 
effect. 

An equally remarkable phenomenon is found in networks with fixed total 
bandwidths per node, where bandwidths per link vanish in the high connectivity limit 
and the relation between the chemical potential and capacity is well defined. For 
sufficiently large total bandwidths, we find a phase transition beyond which the chemical 
potential function has to be described by the Maxwell's construction, implying the 
existence of clusters of balanced nodes having a uniform shortage among them. In the 
case of large but finite connectivity, the Maxwell's construction becomes a linear region 
with nonzero slope, implying a less uniform shortage among the balanced nodes. This 
reflects the benefits of increasing the number of connections in resource allocation. When 
adapted to scaled-free networks, deviations from the Maxwell's construction reveal that 
the presence of hubs is able to homogenize the resources among the nodes with low 
connectivities. 

For future extensions, the theory and algorithms in this paper can be generalized 
to model real networks with inhomogeneous connectivity and bandwidths. Scale-free 
networks with adjustable bandwidth distributions would be one of the most interesting 
systems to study. Further studies can also be done on minimizing the bottleneck effect 
by using heterogeneous bandwidths on different links. It is also worthwhile to consider 
other non-linear shortage costs. We believe that the techniques presented in this paper 
are useful in many different network optimization problems and will lead to a large 
variety of potential applications. 
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